1 概述

转录组 ( Transcriptome ) 是指在特定环境或生理条件下的一个细胞或一个细胞群中所有 RNA 分子的集合，包括信使 RNA ( mRNA )、核糖体 RNA ( rRNA )、转运 RNA ( tRNA )及非编码 RNA ( ncRNA )。它有时被用来指所有的 RNA ，或仅指 mRNA ，指何种 RNA 取决于我们设计的特定实验，以下转录组分析仅针对 mRNA 进行研究。转录组不同于正常基因组的稳定性，它会随外界环境的变化而发生动态变化。因此，转录组可以反映不同实验条件和不同时间里基因表达的变化，但是无法获取转录衰减等 mRNA 降解现象的信息。根据实验设计获得不同的实验样本后，提取 mRNA 并通过二代测序获取 RNA 序列信息用于后续信息分析，本报告分析内容主要以转录组的差异基因分析为重点，并进行功能注释来进一步探究不同的实验设计或环境是如何影响样本的基因表达和功能的发挥

1.1 建库测序

转录组测序一般需要经过 RNA 样品提取，RNA 检测，测序文库构建，上机测序等一系列实验环节，最终得到用于生物信息分析的测序数据。其中每一个环节都会影响测序数据的产出和质量，从而影响后续生物信息分析结果的准确性。镁伽科技重视每一个实验环节的准确性和可靠性，以确保得到高质量的测序数据。具体建库测序流程如下：

建库流程

建库流程

1.2 样品检测

高质量的RNA是整个项目成功的基础，为保证测序数据准确性，我们对 RNA 样品的检测主要包括以下几点：

  • NanoDrop™ One/OneC 检测 RNA 纯度 (OD260/280, OD260/230 比值 )
  • Life Invitrogen Qubit® 3.0 荧光定量仪采用 Qubit™ RNA HS Assay Kit 精确 定量
  • Agilent 4200 TapeStation 系统精确检测 RNA 完整性 ( RIN值 )。 90%，细胞浓度为 700～1200 个细胞/μl

1.3 文库构建与质检

1.3.1 文库构建

样品检测合格之后，进行样品的文库构建，具体mRNA 捕获建库过程如下：

  • 真核生物大部分 mRNA 具有 polyA 的结构，采用带 Oligo ( DT ) 的磁珠对具有 polyA 结构的 mRNA 进行捕获；
  • 以片段化的 mRNA 为模版，随机寡核苷酸为引物，在 M-MuLV 逆转录酶体系中合成 cDNA 第一条链；
  • 用 RNaseH 降解 RNA 链，并在 DNA polymerase I 体系下，以 dNTPs 为原料合成 cDNA 第二条链；
  • 纯化双链 cDNA ，并进行末端修复、加 A 尾和连接测序接头； 用 AMPure XP beads 筛选 200bp 左右的 cDNA，进行 PCR 扩增并再次使用 AMPure XP beads 纯化 PCR 产物，最终获得测序文库。

1.3.2 文库质检

文库构建成功之后对文库进行检测，库检主要包括以下两方面：

  • Kapa qPCR 定量：精确定量文库浓度。
  • Agilent 4200 TapeStation 检测：精确检测文库片段大小。

1.3.3 上机测序

库检合格后，把不同文库按照有效浓度及目标下机数据量的需求 pooling 后进行 Illumina PE150 测序。

2 信息分析流程

目前 mRNA 转录组测序生物信息分析过程主要包括数据质控，参考基因组比对，基因表达量分析，差异基因分析，功能富集分析、差异基因的蛋白互作分析，变异检测、差异可变剪接分析、GSEA分析。分析流程如下图所示：

信息分析流程

信息分析流程

2.1 原始数据质量控制

2.1.1 测序数据说明

测序获得的原始数据 ( Raw Data ) 存在一定比例低质量数据，会对后续信息分析造成很大干扰。高质量的测序数据是后续信息分析结果可靠性的前提，因此需要对原始数据进行严格的数据预处理和评估。我们使用自主研发的开源软件 fastp 1对原始数据进行处理，对低质量序列进行过滤，过滤内容包括如下部分：

  • 去除Reads 中接头 ( Adapter )序列
  • 去除 N ( N 表示模糊碱基信息 ) 碱基数达到一定比例的 Reads (默认 5bp)
  • 去除所含低质量碱基 (质量值 ≤20) 超过一定比例 (默认 40%) 的 Reads
  • 滑窗质量剪裁：以一定大小（默认 4bp ）为滑动窗大小，平均质量低于特定值（默认 20 ）的进行剪裁 上述的处理步骤均是对 Reads 1 和 Reads 2 进行操作。更多详细信息请查看results/qc 文件夹

2.1.2 测序错误率分布

测序过程本身存在机器错误的可能性，测序错误率分布检查可以反映测序数据的质量，序列信息中每个碱基的测序质量值保存在 fastq 文件中。如果测序错误率用 e 表示，Illumina 的碱基质量值用 Qphred 表示，则有：Qphred=-10log10(e) 。Illumina Casava 1.8 版本碱基识别与 Phred 分值之间的简明对应关系见下表。

Phred 分值 不正确的碱基识别 碱基正确识别率 Q.sorce
10 1/10 90% Q10
20 1/100 99% Q20
30 1/1000 99.90% Q30
40 1/10000 99.99% Q40

2.2 参考基因组比对分析

选取STAR （ Spliced Transcripts Alignment to a Reference ）软件2 将转录组测序 Reads 比对到参考基因组。STAR具有速度快、灵敏度和准确度高、内存消耗低等优势。Clean Reads 比对到参考基因组上获得比对特征信息并以 SAM 格式( Sequence Alignment/Map format ) 存储。SAM 序列比对格式标准由sanger制定 3 ，主要用于测序序列比对到基因组上信息特征的展示，BAM 文件是 SAM 文件的二进制格式，存储资源更小，常用于后续计算。基于比对结果使用 Samtools 对 BAM 文件进行排序、过滤重复 Reads 等处理，最后对序列在参考基因组上的比对结果进行统计和绘图展示。

2.2.1 比对结果统计分析

数据下机经过质控处理后，将数据比对到参考基因组上，具体的比对结果统计见文件： results/qc/multiqc_report.html

2.3 表达量分析

表达量分析是差异分析的前提，表达量分析的实质是根据各样本的比对结果统计不同基因的 Reads Count 值，即比对到该基因的 Reads 数有多少，也是后续差异基因分析的基础信息，如果基因的表达量高，则其 Reads 覆盖数也会高。

2.3.0.1 Reads Count 统计表
样本Reads Count统计表
gene N12 N13 N15 T12_D T12_NC T13_D T13_NC T15_D T15_NC
5S_rRNA 6 0 1 0 1 1 5 1 0
5_8S_rRNA 0 0 0 0 0 0 0 0 0
7SK 1 3 2 4 0 3 0 0 0
A1BG 0 1 3 0 0 3 0 1 0
A1BG-AS1 4 2 4 2 1 0 0 2 1
A1CF 1 2 1 1 5 0 0 1 3
  • 表列 ：测序样本 ID
  • 表行 ：各样本中基因 Reads Count 数
  • 所有基因在样本中的 Reads Count 统计请查看 results/counts/all.symbol.tsv

2.3.1 TPM 统计

TPM ( Transcripts Per Kilobase of exon model per Million mapped reads ) 是指每千个碱基的转录每百万映射读取的转录本数，是 RNA-Seq 的双端测序基因表达值的计算方式，先后对基因长度和测序深度进行校正后得到校正值，我们对 Reads Count 进行 TPM 值转换，并绘制 TPM 箱线图，通过图形可以观察不同实验条件下的基因表达水平变化。
TPM箱线图

TPM箱线图

2.3.2 样本相关性

统计学上相关性分析 ( Correlation Analysis ) 是指对两个或多个具备相关性的变量元素进行分析，从而衡量两个变量因素的相关密切程度。对于 RNA-Seq 实验来说，生物学重复的实验样本可以通过相关性分析查看生物学实验操作的变动大小，非生物学重复关系的样本可以通过相关性分析观察不同实验条件间相似性和差异程度。通过各样本的基因定量关系，绘制样本相关性热图，观察各样本间的关系。相关系数越接近 1，说明样品间表达的模式相似度越高。在下图中颜色越深代表相关性越强。
样本TPM表达矩阵相关性热图

样本TPM表达矩阵相关性热图

2.3.3 主成分分析

主成分分析（PCA）是将多个变量通过降维为少数几个相互独立的变量（即主成分），同时尽可能多地保留原始数据信息的一种多元统计分析方法。在转录组的分析中，PCA 将样本所包含的大量基因表达量信息降维为少数几个互相无关的主成分，以进行样本间的比较，图中样本间的距离反应了样本的相似性和离群情况，距离越近，样本越相似。
pca 散点图

pca 散点图

2.4 差异表达基因分析

差异分析 ( Difference analysis ) 是在一组样本中找出产生差异的原因及其对造成差异的影响程度。根据基因的定量分析信息可进行差异基因分析，筛选不同的实验条件下的功能基因。基因差异表达分析的输入数据为基因表达水平分析中得到的Reads Count数据，分析主要分为三部分：

  • 首先对 Reads Count 进行标准化（ normalization ）
  • 然后根据模型进行假设检验的 P-value 以及 FoldChange 的计算；
  • 最后进行多重假设检验校正，得到矫正后的 P 值( padj )
类型 软件 标准化方法 p值矫正方法
有生物学重复 Deseq24 软件自设 BH5
无生物学重复 EdgeR6 软件自设 BH

2.4.1 差异表达基因统计

基因表达具有时间和空间特异性，在两个不同条件下，表达水平存在显著差异，各比较组合间分析得到的差异表达基因统计如下表所示：

差异基因表达列表
gene baseMean log2FoldChange lfcSE pvalue padj
LINC03061 1456.5710438 -2.2579168 0.2199267 0.0000000 0
FSTL1 320.3733667 -3.2547143 0.3367156 0.0000000 0
SCEL 746.2850365 -1.8190108 0.2019606 0.0000000 0
REG1A 495.0132112 5.2736542 0.6889077 0.0000000 0
TMEM238L 294.8064412 1.7754076 0.2359544 0.0000000 0
SERPINB2 186.6521241 -4.4830945 0.5878698 0.0000000 0
DDX11L2 0.3177327 0.0000049 0.0089362 0.8024778 1
WASH7P 128.3431754 0.0001754 0.0089322 0.5105055 1
MIR6859-1 2.5199039 -0.0000302 0.0089360 0.6187630 1
MIR1302-2HG 0.2261730 0.0000000 0.0089362 1.0000000 1
ENSG00000238009 0.3863901 -0.0000056 0.0089362 0.7771276 1
CICP27 3.0549363 -0.0000006 0.0089360 0.9924111 1
  • gene：差异基因的名称
  • baseMean: 所有样本均值
  • log2FoldChange：FoldChange 的 log2 值
  • pvalue：显著性检验的 P值
  • padj ：多重假设检验经过 BH 方法校正后的 P 值

2.4.2 差异表达基因火山图

火山图 ( Volcano Plot ) 是一种展示 FoldChange 和 Pvalue 这两个重要指标的一种统计图，可以非常直观的查看在两样本间发生差异表达的基因整体分布情况。我们根据设定的阈值：padj < 0.05 ，对样本间的差异基因进行筛选获取上调和下调基因相对应的 FoldChange 和 padj 。绘制下图，图中蓝色表示为表达量下调基因，红色表示为表达量上调基因。
差异基因表达火山图

差异基因表达火山图

2.4.3 基因数目韦恩图

韦恩图展示不同组合间差异表达基因的重叠情况，圈内所有数字之和代表该比较组合差异基因总数，重叠区域表示组合间共有的差异基因。对于两个分组(样本)间的比较情况，则展示两组(样本)之间共同表达与特异表达基因的情况，圈内所有数字之和代表该组（样本）表达基因总数，重叠区域表示组(样本)间共有的表达基因数目。
基因表达数目韦恩图

基因表达数目韦恩图

2.5 富集分析

通过富集分析，可以找到不同条件下的差异基因与哪些生物学功能或通路显著性相关，从而揭示和理解生物学过程的基本分子机制。 我们采用 clusterProfiler7 软件对差异基因集进行GO功能富集分析，KEGG通路富集分析等。富集分析基于超几何分布原理，其中差异基因集为差异显著分析所得差异基因并注释到GO或KEGG数据库的基因集，背景基因集为所有参与差异分析并注释到GO或KEGG数据库的基因集。富集分析结果是对每个差异比较组合的所有差异基因集进行富集。

富集分析原理图

富集分析原理图

2.5.1 GO 富集

GO(Gene Ontology)是描述基因功能的综合性数据库（http://www.geneontology.org/)，可分为生物过程（Biological Process）和细胞组成（Cellular Component）分子功能（Molecular Function）三个部分。GO 功能显著性富集分析给出与基因集背景相比，在差异表达基因中显著富集的 GO 功能 条目，从而给出差异表达基因与哪些生物学功能显著相关。该分析首先把所有差异表达基因向 Gene Ontology 数据库的各个 term 映射，计算每个 term 的基因数目，然后找出与整个基因集背景相比，在差异表达基因中显著富集。GO富集分析方法为clusterProfiler的enricher函数，该方法使用超几何分布检验的方法获得显著富集的GO Term。统计被显著富集的各个GOterm中的基因数,GO富集以padj小于0.05为显著富集。

2.5.1.1 柱状图展示
2.5.1.1.1 BP
GO_BP

GO_BP

2.5.1.1.2 MF
GO_MF

GO_MF

2.5.1.1.3 CC
GO_CC

GO_CC

2.5.2 KEGG 富集

在生物体内，不同基因相互协调行使其生物学功能，通过Pathway显著性富集能确定差异表达基因参与的最主要生化代谢途径和信号转导途径。 KEGG(Kyoto Encyclopedia of Genes and Genomes)是有关Pathway的主要公共数据库，其中整合了基因组化学和系统功能信息等众多内容。信息分析时，我们将KEGG数据库按照动物，植物，真菌等进行分类，依据研究物种的种属选择相应的类别进行分析。Pathway显著性富集分析以KEGG Pathway为单位，应用超几何检验，找出差异基因相对于所有有注释的基因显著富集的pathway，KEGG通路富集同样以padj小于0.05作为显著性富集的阈值。

2.6 差异基因蛋白互作网络

蛋白互作网络 ( Protein-Protein Interaction Networks ，PPI ) 是蛋白通过彼此之间的相互作用构成，参与基因表达调节、生物信号传递、细胞周期调控及能量物质代谢等生命过程的各个环节。我们使用 STRING ( https://string-db.org/ ) 蛋白质互作数据库中的互作关系进行差异基因的 PPI 分析，对于差异基因的蛋白互作网路数据文件可使用 StringDB 进行可视化展示8。

2.6.1 差异基因蛋白互作网络结果展示

详细的差异基因蛋白互作网络见文件result/enrichment/ * / * /STRING_info.txt

2.7 差异可变剪接分析

可变剪接 ( Alternative Splicing ) 是指一个 mRNA 前体通过不同的剪接方式 ( 选择不同的剪接位点 ) 产生不同的 mRNA 剪接异构体，是真核生物基因和蛋白质数量有较大差异的重要原因。可变剪接可分为

  • skipped exon (SE)外显子跳跃
  • alternative 5’ splice site (A5SS)第一个外显子可变剪切
  • alternative 3’ splice site (A3SS)最后一个外显子可变剪切
  • mutually exclusive exons (MXE)外显子选择性跳跃
  • retained intron (RI)内含子滞留
等5种类型。我们使用 rmats软件进行样本间的差异可变剪接分析，可变剪接形式和软件检测流程原理如下图所示，主要根据 altered junction reads 和 altered coverage of AS-exons 计算可变剪接的差异显著性。
可变剪接分析原理

可变剪接分析原理

2.7.1 差异可变剪接统计

每个比较组合得到的差异可变剪接统计如下表：

差异可变剪接统计表
EventType EventTypeDescription TotalEventsJC TotalEventsJCEC SignificantEventsJC SigEventsJCSample1HigherInclusion SigEventsJCSample2HigherInclusion SignificantEventsJCEC SigEventsJCECSample1HigherInclusion SigEventsJCECSample2HigherInclusion
SE skipped exon 83845 83845 0 0 0 0 0 0
A5SS alternative 5’ splice sites 9546 9546 0 0 0 0 0 0
A3SS alternative 3’ splice sites 14512 14512 0 0 0 0 0 0
MXE mutually exclusive exons 8716 8716 0 0 0 0 0 0
RI retained intron 9881 9881 0 0 0 0 0 0
SE skipped exon 85040 85040 0 0 0 0 0 0
A5SS alternative 5’ splice sites 9613 9613 0 0 0 0 0 0
A3SS alternative 3’ splice sites 14564 14564 0 0 0 0 0 0
MXE mutually exclusive exons 8975 8975 0 0 0 0 0 0
RI retained intron 9895 9895 0 0 0 0 0 0
  • EventType：可变剪接类型;
  • EventTypeDescription：可变剪接类型描述;
  • TotalEventsJC ：跨越剪切点的事件总数;
  • TotalEventsJCEC：跨越剪切点但不跨越外显子边界的事件总数;
  • SignificantEventsJC：显著的JC总数;
  • SigEventsJCSample1HigherInclusion：实验组显著的JC数量;
  • SigEventsJCSample2HigherInclusion：对照组显著的JC数量;
  • SignificantEventsJCEC：显著的JCEC总数;
  • SigEventsJCECSample1HigherInclusion：实验组显著的JCEC数量;
  • SigEventsJCECSample2HigherInclusion：对照组显著的JCEC数量。

2.8 GSEA分析

GSEA分析（Gene Set Enrichment Analysis），基因集富集分析9，是目前非常常用的RNA表达分析手段。GSEA不需要指定明确的差异基因阈值，算法会根据实际数据的整体趋势，为我们提供一种合理的解决方法，即使在没有先验经验的情况下也能在表达谱整体层次上对多个基因进行分析，从而从数理统计层面把表达谱数据与生物学意义很好地衔接起来。目前分析物种只支持人、小鼠，最低条件满足3V3组合才能跑gsea分析。

2.8.1 GSEA结果展示

下图展示了基因表达的可变性与细胞中平均表达水平的关系，红线显示的了基于此关系的期望值， 我们将用于后续的细胞亚群识别或沿轨迹排序细胞标记为黑点，而其他基于显示为灰点。

GSEA结果展示表
X ID Description setSize enrichmentScore NES pvalue p.adjust qvalue rank leading_edge core_enrichment
hsa04390 hsa04390 Hippo signaling pathway 129 -0.6189196 -2.311827 0.0000004 0.0000234 0.0000107 1187 tags=33%, list=8%, signal=31% 7529/6788/7533/659/5518/10971/3398/23291/324/4771/4087/126374/7004/5515/7046/7532/2932/5501/896/5516/25937/26524/658/55233/7534/999/60/332/7477/8323/154796/84962/71/2736/1490/4092/1739/81029/23286/7476/595/894/3689
hsa03010 hsa03010 Ribosome 129 -0.6128998 -2.289342 0.0000008 0.0000234 0.0000107 3903 tags=70%, list=26%, signal=52% 63931/29074/6231/11224/9349/6147/6135/9045/6139/6132/6222/6156/6154/6169/219927/6157/6168/6165/29093/6142/7311/6188/6209/6183/6208/51121/6217/51264/6143/6161/6227/6155/6160/6166/51263/6228/6130/6230/6171/11222/6133/6152/6146/6234/6235/6137/6176/6129/6189/6218/55052/6122/51187/6136/6187/51318/6125/6210/6181/6150/6204/6134/51021/6138/6124/6201/6158/28998/6175/3921/64928/6232/6223/6164/6203/4736/29088/55173/6224/9553/23521/9801/200916/6193/64965/6144/6207/388524/6233/6206
hsa04110 hsa04110 Cell cycle 154 -0.4910432 -1.864148 0.0005038 0.0104127 0.0047732 4570 tags=57%, list=30%, signal=40% 9184/5519/10926/991/10735/23595/151648/4176/27085/1387/23047/113130/4089/51434/25/81620/3065/993/4998/4171/8317/9978/9088/890/8318/9134/1019/8697/4174/57082/8454/4175/4173/51529/1031/701/2810/5347/1871/983/891/9133/10393/5885/64682/7531/1870/6502/3066/7029/4085/1869/5001/5111/8379/6500/5591/23244/9319/10403/5925/5527/7529/545/7533/23063/5518/990/10971/1111/5529/157570/4087/5515/7532/25847/5933/2932/896/5516/7027/7534/7157/10274/1021/1029/595/894
hsa05132 hsa05132 Salmonella infection 213 -0.4463716 -1.772945 0.0006940 0.0107566 0.0049309 3402 tags=40%, list=22%, signal=32% 6416/9367/1499/1432/10376/1781/5291/26999/6188/4074/1434/203068/4790/5605/1639/5868/302/5216/8767/10552/5585/5606/51143/5878/10392/140735/207/5287/836/10121/65082/6885/3840/10787/84516/3265/3836/5898/375/10093/5601/10006/6500/81873/382/6281/27072/55207/1778/84790/3799/55770/5594/3839/8976/10096/998/387/8772/1147/5420/11258/10671/9842/10097/10640/112574/54205/4615/10627/8655/103910/55860/840/3654/60/55845/79026/5288/10092/71/9265/7846/3606/2316/10398
hsa04360 hsa04360 Axon guidance 144 -0.4904220 -1.846516 0.0012931 0.0111653 0.0051182 1395 tags=24%, list=9%, signal=22% 3688/5534/5594/818/4233/23365/659/8440/5361/998/387/5295/1072/2773/23654/3845/5062/10627/1946/2932/1948/658/103910/4775/29984/4893/5781/1969/1808/10371/81029/8503/56963/10512/10398
hsa04530 hsa04530 Tight junction 142 -0.4852686 -1.826337 0.0014407 0.0111653 0.0051182 1689 tags=30%, list=11%, signal=27% 5111/10093/5601/81873/7074/3996/3688/84790/1366/10096/998/5518/387/153562/81/11346/10097/4771/5515/10627/5516/1364/100506658/7082/103910/60/154810/137075/79784/51762/10092/154796/4627/71/9076/7846/7430/1739/51421/5581/595/10398
hsa04510 hsa04510 Focal adhesion 186 -0.5988323 -2.061757 0.0000000 0.0000000 0.0000000 3027 tags=31%, list=12%, signal=27% 3918/5595/4659/3655/5291/8394/1398/71/7408/29780/10627/6655/6093/3611/3694/10398/5747/5500/3673/3678/2013/3696/5604/4233/5501/5170/8395/673/5829/56034/2909/3691/5599/5881/3675/2317/3909/3913/5906/998/3688/7414/6654/329/2932/331/4660/595/824/2316/1956/5601/2064/340156/5908/394/7791
hsa04810 hsa04810 Regulation of actin cytoskeleton 206 -0.5387199 -1.869518 0.0000000 0.0000006 0.0000003 3757 tags=34%, list=15%, signal=29% 4628/1793/5216/26999/27006/9459/10000/8826/2254/10097/5962/10109/5595/4659/3655/5291/85464/8394/1398/71/10152/10096/10627/6655/6093/3694/54434/3984/10092/10398/3845/5747/10163/5500/3673/3678/2934/3696/5604/8936/55740/6548/28964/5501/8395/673/5829/56034/1730/2909/3691/81624/5881/3675/324/998/3688/7414/6654/10297/4660/55970/4629/1956/23365/55845/200576/340156/79784/4627/85477
hsa04820 hsa04820 Cytoskeleton in muscle cells 203 -0.5145934 -1.786836 0.0000008 0.0000230 0.0000118 4996 tags=43%, list=20%, signal=35% 25802/58/482/2026/57731/3685/2274/6443/287/7139/1158/6385/7170/29765/6382/23224/3674/91977/89/83660/3908/3695/93661/2010/64236/7273/9499/9672/481/1829/23002/1292/3672/6712/7111/1277/57644/2192/4628/832/1605/29767/9260/3339/6641/1824/8557/8910/79188/3728/3655/6709/22989/478/5318/9124/6711/71/7168/476/271/161176/3694/10611/10398/3673/3678/3696/1466/1730/829/3691/81624/3675/3688/7414/7134/288/4629/6645/29766/2275/7171/1832/7791/79784/4627
hsa04530 hsa04530 Tight junction 158 -0.5167376 -1.742798 0.0000097 0.0002095 0.0001077 3782 tags=39%, list=15%, signal=34% 5566/4628/4217/57826/10097/51422/5962/5520/9080/10109/55114/4771/9076/861/9414/2017/9693/51735/5581/51421/71/3993/7408/10207/10096/9368/10627/6093/10092/10398/54566/9223/50848/5562/11346/9181/5563/55844/100506658/84952/5584/57530/56288/5599/92359/123720/5906/154796/998/3688/4301/137075/595/1739/4629/5601/2064/4218/79784/5590/4627/51762
hsa05165 hsa05165 Human papillomavirus infection 301 -0.4410204 -1.587632 0.0000145 0.0002497 0.0001284 3512 tags=33%, list=14%, signal=29% 4855/8324/3716/1277/6772/245973/1027/3516/9134/5566/4854/8325/10000/6199/3455/5734/83439/1387/84441/1385/55012/8717/894/6009/3280/148022/5520/3918/10379/1857/5315/5595/182/4600/7098/3655/64764/5291/7976/472/5528/535/6773/528/5663/3993/54626/10207/6198/1017/9368/1021/6655/51382/5527/3694/6934/3845/5747/9223/3673/3678/3696/5604/5610/23352/55844/3066/7132/5829/9794/29110/5584/3691/56288/7476/5529/92359/5525/3675/3909/2033/3913/324/998/3688/8638/3065/8772/6654/10297/2932/595/3454/1739/1956/5523/7855/5590/55534
hsa04144 hsa04144 Endocytosis 243 -0.4614140 -1.630226 0.0000177 0.0002538 0.0001305 3471 tags=33%, list=14%, signal=28% 1211/1785/11031/10097/3303/9897/22841/60682/10109/3482/58513/9525/2060/50807/55040/23624/4087/116985/10564/378/10193/9765/8411/7048/8394/8724/84440/10096/137492/11267/387680/55048/10015/10092/57403/10617/382/7879/3304/10565/28964/11059/55616/51534/116987/8395/23325/57132/30844/5584/79720/829/8853/56288/9372/9135/408/9815/3949/6642/55680/998/80223/26286/80230/64744/83737/9265/5868/9267/116984/1956/30846/4218/64145/5869/5590/8766/9101
  • ONTOLOGY：GO类别
  • ID：GO编号
  • Description：GO描述
  • setSize：该通路中包含表达数据集的基因数目
  • enrichmentScore：富集分数
  • NES：标准化后的富集分数
  • pvalue ：统计学差异显著性检验的P值，一般为 0.01 或 0.05
  • p.adjust ：多重假设检验经过 BH 方法校正后的 P 值
  • qvalue: Q检验法矫正后的p值
  • rank：当 ES 最大时，对应基因所在排序好的基因列表中所处的位置
  • leading_edge：tags 表示核心基因占该通路基因集的百分比；list 表示核心基因占所有基因的百分比；signal，将前 2 项统计值结合在一起计算出的富集信号强度

2.9 免疫浸润分析

组织是由不同谱系和亚型的细胞类型组成的复杂环境，每种细胞都有自己独特的转录组。因此，批量转录组分析是细胞类型特异性基因表达的总和加权的细胞类型比例在给定的样本。去卷积的基因表达谱允许重建组织的细胞组成。xCell 是一个强大的计算方法，转换基因表达谱为丰富分数的免疫和基质细胞类型跨样本。

2.9.1 xCell浸润分析结果

免疫浸润分数热图

免疫浸润分数热图

3 软件列表

分析内容 软件 版本
差异分析 deseq2 1.38
差异分析 edgeR 3.40.2
富集分析 clusterProfiler 4.10.0
剪接分析 rmats 4.3.0
数据质控 rseqc 5.0
免疫浸润分析 xcell 1.3

4 参考文献

  1. Shifu, C. , Yanqing, Z. , Yaru, C. , & Jia, G. (2018). Fastp: an ultra-fast all-in-one fastq preprocessor. Bioinformatics, 34(17), i884-i890.↩︎

  2. Dobin, A., et al. (2013). STAR: ultrafast universal RNA-seq aligner. Bioinformatics (Oxford, England), 29(1), 15–21. https://doi.org/10.1093/bioinformatics/bts635↩︎

  3. Li, H. E. A. (2009). The sequence alignment / map ( sam ) format. Bioinformatics, 25(1 Pt 2), 1653-4.↩︎

  4. Love, M. I. , Huber, W., & Anders, S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol, 15(12), 550.↩︎

  5. Benjamini, Y. , Hochberg, Y. (1995). Controlling the False Discovery Rate - a Practical and Powerful Approach to Multiple Testing. J Roy Stat Soc B, 57(1), 289-300.↩︎

  6. Robinson, MD. , McCarthy, DJ. , Smyth, GK. (2010). edgeR: a bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics, 26, 139–140.↩︎

  7. Yu, G. , Wang, L. G. , Han, Y. , & He, Q. Y. (2012). Clusterprofiler: an r package for comparing biological themes among gene clusters. Omics-a Journal of Integrative Biology, 16(5), 284-287.↩︎

  8. Szklarczyk D, Gable AL, Nastou KC, Lyon D, Kirsch R, Pyysalo S, Doncheva NT, Legeay M, Fang T, Bork P, Jensen LJ, von Mering C (2021). “The STRING database in 2021: customizable protein-protein networks, and functional characterization of user-uploaded gene/measurement sets.” Nucleic Acids Research (Database issue), 49.↩︎

  9. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A. 2005;102:15545–15550.↩︎